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Stability  Analysis  for 
Rotating  Stall  Dynamics 
in  Axial  Flow  Compressors* 
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Abstract.  This  paper  analyzes  the  nonlinear  phenomenon  of  rotating  stall  via  the  methods 
of  projection  [Elementary  Stability  and  Bifurcation  Theory,  G.  Iooss  and  D.  D.  Joseph, 
Springer- Verlag,  1980]  and  Lyapunov  [J.-H.  Fu,  Math.  Control  Signal  Systems,  7,  255- 
278,  1994],  A  compressor  model  of  Moore  and  Greitzer  is  adopted  in  which  rotating  stall 
dynamics  are  associated  with  Hopf  bifurcations.  Local  stability  for  each  pair  of  the  critical 
modes  is  studied  and  characterized.  It  is  shown  that  local  stability  of  individual  pairs  of 
the  critical  modes  determines  collectively  local  stability  of  the  compressor  model.  Explicit 
conditions  are  obtained  for  local  stability  of  rotating  stall  which  offer  new  insight  into  the 
design,  and  active  control  of  axial  flow  compressors. 


1.  Introduction 

Axial  flow  compressors  are  vital  parts  of  the  aeroengines.  Yet  two  distinct  aerody¬ 
namic  instabilities,  rotating  stall  and  surge,  can  severely  limit  compressor  perfor¬ 
mance.  Both  of  these  instabilities  are  disruptions  of  the  normal  operating  condition, 
which  is  designed  for  steady  and  axisymmetric  flow.  Rotating  stall  is  a  severely 
nonaxisymmetric  distribution  of  axial  flow  velocity  around  the  annulus  of  the 
compressor,  taking  the  form  of  a  wave  or  stall  cell  that  propagates  steadily  in  the 
circumferential  direction.  Surge,  on  the  other  hand,  is  an  axisymmetric  oscillation 
of  the  mass  flow  along  the  axial  length  of  the  compressor.  Both  of  these  disruptions 
can  have  catastrophic  consequences  for  jet  airplanes. 
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The  nonlinear  phenomena  of  rotating  stall  and  surge  have  been  studied  exten¬ 
sively  in  the  literature  in  which  a  PDE  (partial  differential  equation)  model  of 
Moore  and  Greitzer  [  10]  becomes  standard.  In  this  PDE  model,  rotating  stall  dy¬ 
namics  are  described  by  a  PDE,  and  surge  is  described  by  two  ODEs  (ordinary 
differential  equations).  If  only  the  fundamental  harmonics  in  a  spatial  Fourier  se¬ 
ries  of  the  disturbance  flow  is  retained,  the  PDE  model  reduces  to  a  third-order 
ODE  model.  The  bifurcation  phenomena  associated  with  rotating  stall  and  surge 
for  this  model  are  analyzed  thoroughly  in  [2],  [9]. 

Because  the  approximation  error  is  due  to  the  truncation  of  high-order  harmon¬ 
ics  in  the  disturbance  flow,  this  paper  will  focus  on  stability  analysis  for  rotating 
stall  dynamics.  The  original  PDE  model  of  [10]  will  be  investigated  indirectly 
through  transformation  into  an  ODE  model  of  order  2 N  +  2,  where  N  -*  oo  is 
the  number  of  pairs  of  stall  modes.  This  multimode  model  ( N  >  1)  gives  much 
more  accurate  representation  and  recovers  the  PDE  model  when  N  ->  oo.  The 
equilibria  of  the  uniform  flow  will  be  derived,  and  the  rotating  stall  dynamics  will 
be  characterized  for  the  multimode  model.  The  local  stability  of  the  compressor 
model  will  be  studied  via  both  the  methods  of  projection  [7]  and  Lyapunov  [3], 
for  which  necessary  and  sufficient  conditions  are  obtained  in  terms  of  various 
compressor  parameters. 

Our  results  show  that  the  stability  condition  for  each  pair  of  the  critical  modes 
determines  the  stability  for  rotating  stall  collectively.  We  point  out  that  although 
in  a  special  case  our  stability  condition  for  rotating  stall  dynamics  is  the  same  as 
that  for  the  simple  third-order  ODE  model,  stability  analysis  for  the  third-order 
ODE  model  cannot  replace  that  for  the  PDE  or  for  the  multimode  model.  In  fact, 
active  control  for  rotating  stall  based  on  the  third-order  ODE  model  often  fails  to 
stabilize  the  critical  operating  point  of  the  PDE  model  because  of  the  neglected 
high-order  harmonics  in  the  disturbance  flow. 

Interested  readers  are  referred  to  [4]  for  a  more  complete  reference  on  stability 
and  stabilization  of  the  rotating  stall  dynamics  for  the  third-order  Moore-Greitzer 
model.  Our  stability  analysis  for  the  multimode  model  with  N  — >  oo  and  the 
explicit  conditions  for  local  stability  of  rotating  stall  offer  new  insight  into  the 
design,  and  active  control  for  axial  flow  compressors. 

A  schematic  axial  flow  compressor  is  shown  in  Figure  1 .  The  following  notation 
is  used  in  this  paper: 


R  = 

mean  rotor  radius 

lc  = 

h  +  Is  +  1/a 

Ac  = 

compressor  duct  area 

vP  = 

volume  of  plenum 

U= 

blade  speed  at  mean  radius 

IM  = 

viscosity  coefficient 

as  = 

speed  of  sound 

Ps  = 

static  pressure  in  plenum 

B  = 

(U /2as)y/vpf{AcLc) 

Pr  = 

total  pressure  ahead  of  entrance 

W  = 

semiwidth  of  cubic  characteristic 

and  following  the  throttle  duct 

a,b  = 

time  lag  of  blade  passage 

*  = 

(ps  —  pr  )jpU 2 :  pressure  rise 

Rotating  Stall  Dynamics  333 


H  = 

semiheight  of  cubic  characteristic 

<P  - 

local  flow  coefficient  at  station  0 

Y  = 

throttle  position 

<j>  = 

mean  flow  coefficient  at  station  0 

m  = 

exit  duct  length  factor 

V  = 

disturbance  flow  at  r/ 

h,  h,  h  = 

length  of  exit,  entrance,  throttle 

V  = 

axial  distance  from  station  0 

ducts,  in  wheel  radius 

?  = 

Ut/R  where  t  is  time 

Figure  1.  Schematic  of  compressor  showing  nondimensionalized  lengths. 


2.  Multimode  Moore-Greitzer  model 

A  post-stall  model  for  compression  systems  was  developed  by  Moore  and  Greitzer. 
Its  full  PDE  form  is  described  by  [10] 


where  cp  =  <P  +  <pjf7=0  is  the  local  mass  flow  at  station  0.  The  incompressible 
and  irotational  assumption  on  the  gas  flow  implies  the  existence  of  a  disturbance 
flow  potential  that  satisfies  Laplace’s  equation  with  zero  boundary  condition  at 
=  —If  <  —lc,  which  in  turn  implies  that  the  disturbance  flow  at  station  0  has 
the  form: 

N 

0  =  ( An  cos(n6))  +  Bn  Sin(«0)) ,  N  ->■  oo.  (3) 

n= 1 

Let  {ei0n  }^0  be  a  set  of  2 N  +  1  uniformly  distributed  samples  on  the  unit 
circle.  Denote 

=  <P(02)  HOin+i)], 


(4) 
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i 

V2 

cos  6i 

l 

COS  02 

-1-  i 

■Jl  1 

COS02N+1 

T~  1  2 

sin@i 

sin@2 

sin02A?+i 

V  2N  +  1 

cos  N6\ 

cos  N62  , . . 

• 

cos  N02n+i 

-  sin  N6\ 

sin  NQ2  . . . 

sin  N92N+i  - 

E  =  diag  (Zc,  m\h,  m2I2, . 

. . ,  m^h) . 

m  coahOilf  )  1 

mn  =  — ■  "V  , . t  r  +  -■ 

n  smh(n//r)  a 


h  = 


F  ■=  diag  0, 


l 

L  b 


-a  -i"1 
f1  b 

“Mj 


1  0 
0  1 

-fxN2 

N 

b 


N  -l 
b 

~HN2_ 


)■ 


M\  =  TrE~l  FT  M2  =  TtE~1T. 


(5) 

(6) 
(7) 


(8) 

(9) 


Then  the  full  PDE  model  can  be  approximated  by  the  following  multimode  Moore- 
Greitzer  model  [8]: 


(p n 

4>N  =  M\<Pn  +  M2^c(A<Pn)  -  M2e„4},  &<Pn  =  -  eN 


w 


V  = 


1 


4B2lc  \  2N  -f  1 


yM. 


-T 


=  [1  1 


1]. 


(10) 

(11) 


This  multimode  model  converges  to  the  PDE  model  described  in  (1)  and  (2)  uni¬ 
formly  as  N  oo  under  some  rather  mild  assumptions.  Note  that  =  T{  , 
and 

TeN  =  eNVW+l,  TTeN  =  Jl—.,  (12) 

V2  N  +  1 

and  ej/  =  [  1  0  •  ■  ■  0  ] .  The  performance  characteristic  curve  is  assumed 

to  be  a  cubic  polynomial  [10]  and  is  given  by 

fc(M  =  H  \coeN  +  ci4>n  +  cz<t>N J  ,  (13) 

where  x'2  =  jr'*(fc=2  —  x  ■  x  denotes  a  Hadamard  product  [5],  where  x'k  raises 
each  element  of  x  to  the  power  of  k.  The  uniform  flow  can  be  represented  by 
tpff  =  <pee^,  where  <pe  is  the  intensity  of  the  flow  rate. 

Suppose  that  the  flow  is  perturbed  near  the  equilibria  of  the  uniform  flow.  Then 

4>N  =  +  8<Pn,  &<Pn  =  tyS(pN  +  ^ ~  -  1^  ZN,  (14) 

where  </>e  ^  0  is  a  scalar  factor,  representing  the  flow  rate  intensity.  It  follows  that 


&<Pn  — 


^Mi  +  M2H 


ci  3  c3 
IV  +  W 


ScpN 


3Hc3 

W2 


(|-1  + 
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+ 


(m^  +  M2H 


Co  +  Cl 


epj  —  M2eN}i>- 
(15) 


Proposition  2.1.  The  steady  equilibria  (4>n  ,  'R)  =  (<peen,  xPe)for  the  multimode 
compressor  model  of(  10)  and  (11)  are  dependent  on  the  throttle  position  parameter 
y  and  are  governed  by 


V,  =  H 


Co  +  Cl 


(M 


+  C3 


W 


i/2 


(16) 


Proof.  Consider  the  perturbed  multimode  model  for  <p^  =  <peeN  +  <5 <Pn-  Then 
8(f> at  =  0  is  an  equilibrium  if  and  only  if  there  exists  a  real  such  that 

{^M\<t>e  +  M2H 

(17) 

Denote  4  as  the  identity  matrix  of  size  k.  Substituting  the  expressions  of  M\  and 
Mi  gives 


CO  +  C!  (£-l)+<*(£-l) 


en  —  M2e^^e  =  0. 


TrE~l 


(18) 


Because  the  ( Disposition  of  F  is  zero,  this  equation  is  equivalent  to  the  first 
equality  in  (16).  By  setting  'i'  =  0,  the  equilibrium  point  has  to  satisfy  an 
additional  condition, 


-T 

Ys/^e  =  Six  ,  =  <Pe,  (19) 

2N  +  1 

which  is  equivalent  to  the  second  equality  of  (16).  Hence  (16)  characterizes  the 
equilibria  of  the  compressor  model.  □ 


Proposition  2. 1  indicates  that  the  equilibria  of  the  multimode  model  are  func¬ 
tions  of  y,  the  throttle  parameter,  which  can  in  turn  be  determined  by  the  equilibria 
as  in  (19).  Equation  (19)  is  also  referred  to  as  the  throttle  curve  as  shown  in  Figure  2 
(the  dashed  lines  A-B  and  C-D  for  two  different  y  values). 

Denote  the  state  variable  as  x  —  &<j>N  ®  SW  where  ®  denotes  the  direct  sum 
and  =  ip  —  A  local  model  around  the  equilibria  ((peeN,  'Pe)  has  the  form 


x  —  f  (y,  x)  =  Lx  +  Q[x,  x]  +  C[x,  x,x]  -| - , 


(20) 
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where 


tee 


N 


-tceN  E-V2(F+aI)E-V2 


S, 


Qlx, 


,X]=\  ^  (w  ~  l)  M2  02W+1 


C[x,  x,  x]  — 

In  the  preceding  expressions 
S 


0T 

u2N+l 


x 

4 <l>. 


[^2  02N+1  1  .3 

r  ■*  ■ 

J 


2W+1 


8>KZ 


0  2B-J(2N  +  l)/c 

jBl/2r  o 

i  y 


K  = 


2Blc ' 

a  =  Hir'c  =  H 


c\  +  3c3 


w 


- 1 


8fl2W^’ 

,dfh  =  H_ 
d<t>e  W 

By  the  form  of  e/v,  the  linear  matrix  can  be  written  as 

L  =  S~l DS,  D  =  diag  (D0,  D\,  ■  ■  ■ ,  Dn)  , 

r 

Do  = 


'  r  k 

,  D„  = 

a  —  ;u,n2 

—n/b 

_  — ic  a/lc_ 

n/b 

a  -  p.n2 

(21) 

(22) 

(23) 

(24) 

(25) 

(26) 

(27) 

(28) 


for  1  <  n  <  N.  Because  5  is  a  matrix  of  similarity  transformation,  an  eigenvalue 
X  of  L  satisfies 


det 


X  —  r 


-ice 


N 


iceN  XI  -  E~l/2(F  +  aI)E-V2 


=  0. 


Proposition  2.2.  Consider  the  compressor  model  of  (10)  and  (11).  Its  linearized 
system  near  (<j>eetf,  ^e)  as  in  ( 16)  has  (N  +  1)  pairs  of  eigenvalues  that  are  given 
as  follows: 

A-1,2  -  0.5(r  +  a/lc)  ±  0.5;^4(k2  +  ra/lc)  -  (r  +  a/lc)2,  (29) 
X2n+i,2n+2  =  (a  -  tin2  ±  j^j  m~l ,  (30) 

where  n  =  1 , 2, . . . ,  N,  and  j  —  a/— 1. 


Proof.  Careful  observation  of  SLS  1  with  L  and  S  as  given  in  (21)  yields 


D  =  SLS  1  =diag(D0,Di,...,Z)jv),  D0 
" 2  —n/b 


x  K 
— k  a/lc 


Dn  = 


a  —  jin 

n/b  <x  —  /in 2 


j  mn  l,  n  —  1,  2, 


■,N. 


(31) 

(32) 


Rotating  Stall  Dynamics  337 


Direct  calculation  for  the  eigenvalues  of  L  gives 


N 

det(A7  —  L)  =  det(A7  -  D)  =  f]  det(A7  ~  A.)  =  0, 

n— 0 


<==>•  det(A7  —  Do)  =  A2  —  (r  +  a/lc)  A.  +  k2  +  r  a/lc  =  0, 


det(A/  —  Dn )  =  (m„A  -  a  +  fin2)2  +  n2 /b2  =  0. 

Hence  the  expressions  for  2 (N  +  1)  eigenvalues  of  L  are  those  shown  in  (29) 
and  (30).  □ 


Figure  2.  Schematic  compressor  characteristic,  showing  rotating  stall. 

Let  us  first  examine  the  performance  curve  i/Ac(-)  for  the  special  case  of  ji  —  0. 
See  the  A-D  curve  in  Figure  2,  where  the  dashed  lines  are  the  throttle  curve 
‘Pe  =  d>2/y2,  whose  intersection  with  ijAc(-)  determines  the  operating  point.  The 
maximum  pressure  rise  takes  place  at  a  =  0  (point  A),  where  the  derivative  of 
equals  zero.  Because  a  <  0  on  the  right  side  of  the  maximum  pressure  rise, 
the  N  pairs  of  eigenvalues  of  the  L  matrix  as  given  in  (30)  are  stable.  By  r  <  0 
the  pair  of  eigenvalues  of  (29)  are  also  stable.  Hence  the  compressor  model  is 
locally  stable  at  the  right  of  point  A,  and  the  axisymmetric  flow  <Pn  =  <pe^N  is  a 
stable  equilibrium.  However,  if  the  throttle  value  y  decreases,  then  the  flow  rate 
intensity  <j>e  decreases,  and  the  derivative  of  eventually  becomes  positive. 
Thus  the  compressor  model  is  unstable  at  the  left  of  point  A.  In  this  case  the  N 
pairs  of  eigenvalues  in  (30)  cross  the  imaginary  axis  simultaneously.  This  is  called 
the  critical  case,  for  which  Hopf  bifurcations  occur  that  induce  rotating  stall.  For 
N  =  1,  the  compressor  operates  along  the  stall  curve  rjrs{ •),  as  shown  in  the 
B-C  branch  of  Figure  2  which  is  locally  stable.  If  the  underlying  bifurcations  are 
subcritical,  or  unstable,  the  A-C  portion  of  the  stall  curve  is  unstable.  Stall  cells 
will  be  bom  at  point  A  and  will  grow,  which  will  throttle  the  operating  point  from 
A  to  B  quickly.  There  is  a  tremendous  drop  in  both  the  pressure  rise  and  the  flow 
rate.  Moreover,  increasing  the  throttle  position  and  flow  rate  at  point  B  does  not 


338  Gu,  Sparks,  and  Belta 


increase  the  pressure  rise.  Rather,  the  pressure  rise  becomes  even  Sower  before  it 
reaches  point  C,  at  which  it  jumps  back  to  (again  because  of  loss  of  stability)  the 
performance  curve  i^c(-).  The  hysteresis  loop  A-B-C-D  is  the  main  cause  of  loss  of 
compressor  performance  and  potential  damage  to  aeroengines.  On  the  other  hand, 
if  the  underlying  bifurcations  are  stable,  or  supercritical,  then  point  A  represents 
a  stable  equilibrium.  Decreasing  the  flow  rate  will  not  throttle  the  operating  point 
to  B.  For  N  =  1,  the  pressure  rise  will  gradually  decrease  along  the  curve  A-C 
(which  will  be  above  the  throttle  line  at  y  —  yc)  and  can  thus  be  easily  reversed  by 
increasing  the  throttle  position.  Hence  the  stability  of  the  Hopf  bifurcation  plays 
an  important  role  in  analyzing  the  stability  of  the  rotating  stall  dynamics. 

For  the  case  p.  >  0,  the  critical  value  of  the  flow  rate  intensity  (f>e  takes  place 
on  the  left  of  the  maximum  pressure  rise  because  the  viscosity  coefficient  p  tends 
to  damp  out  the  flow  instability,  thus  delaying  the  occurrence  of  rotating  stall. 
A  further  decrease  in  the  flow  rate  intensity  will  also  cause  instability  of  the 
pair  of  eigenvalues  in  (29)  that  corresponds  to  surge  dynamics.  The  next  result 
characterizes  the  critical  throttle  parameter  corresponding  to  rotating  stall. 


Corollary  2.3.  Suppose  that  c\  +  3c3  =  0  for  the  performance  characteristic 
curve  and  ii  =  0.  Then  the  N  pairs  of  eigenvalues  as  in  (30)  cross  the 
imaginary  axis  at  the  same  critical  throttle  parameter  value 


Y  —  Yc  — 


4>e 


4>e=2W 


2  w  _  2  W 

ff\/rc(2.W)  Jti(co  +  ci+xY 


(33) 


and  A.2n+l,2n+2  =  ±  for  1  <  n  <  N.  For  p  >  0,  the  N  pairs  of  eigenvalues 
as  in  (30)  cross  the  imaginary  axis  at  N  different  critical  throttle  parameter  values 
that  are  given  by 


Y  -  Yen 


Hco  +  HciJl  +  lgjg  +  Hci 


(34) 


and  expression  X 2n+i,2n+2  =  ±  jj,mn  1  remains  the  same  for  n  —  1, 2, . . . ,  N. 


Proof.  Because  ci  +  3c3  =  0,  <t>e  =  2W  is  the  critical  flow  rate  intensity  if  p  —  0. 
Substituting  into  (16)  and  (19)  yields  the  critical  throttle  parameter  value  as  in  (33). 
For  p  jh  0,  the  nth  pair  of  eigenvalues  cross  the  imaginary  axis  at  a  =  fin2,  or 


(35) 


Substituting  into  (16)  yields  the  critical  throttle  parameter  value  as  in  (34).  □ 
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For  /i  >  0,  the  largest  N  value  such  that  yCn  is  real  for  n  =  1,2,...,  N  is 
bounded  as  N  <  -J—3c'sH/{piW).  This  is  the  same  as  in  [2], 3  For  /x  =  0,  the 
largest  N  value  can  be  unbounded.  Moreover  the  critical  values  of  <pe  and  fiL  under 
the  condition  of  ci  +  3c3  =  0  are  given  by 

<pc  =  2W,  *c  =  H(c0  +  ci  +  ca).  (36) 


3.  Stability  analysis  for  rotating  stall 

This  section  considers  rotating  stall  dynamics  with  N  pairs  of  critical  modes  cor¬ 
responding  to  N  harmonics  of  spatial  Fourier  series  for  the  local  flow  rate.  Local 
stability  for  each  pair  of  the  critical  modes  will  be  analyzed  first  using  the  projec¬ 
tion  method  [7],  [6]  (as  outlined  in  the  Appendix).  It  is  known  that  even  if  every 
pair  of  critical  modes  corresponding  to  the  rotating  stall  is  stable,  the  rotating  stall 
dynamics  still  may  not  be  stable  because  of  the  coupling  of  different  pairs  of  crit¬ 
ical  modes.  Hence  the  Lyapunov  method  of  Fu  [3]  (as  outlined  in  the  Appendix) 
will  be  adopted  to  analyze  further  the  local  stability  of  the  rotating  stall.  It  will  be 
demonstrated  that  for  the  case  of  /x  =  0,  local  stability  of  each  pair  of  the  critical 
modes  implies  local  stability  of  the  rotating  stall  dynamics. 

From  Section  2,  a  local  model  of  the  multimode  Moore-Greitzer  model  has  the 
form  in  (20),  where  linear,  quadratic,  and  cubic  terms  are  given  as  in  (21)  through 
(23).  In  order  to  apply  the  projection  and  the  Lyapunov  method,  left  and  right 
critical  eigenvectors  corresponding  to  the  rotating  stall  must  be  computed.  Denote 
©„  as  a  column  vector  of  size  2 N  +  I  as 

en  =  [e~Jn91  e~inei  ■■■  e~in9w+'  ]T  ,  (37) 

where  the  ejdk  ’s  are  equally  distributed  on  the  unit  circle  as  assumed. 


Lemma  3.1.  The  right  and  left  eigenvectors  corresponding  to  the  nth  pair  of  the 
critical  eigenvalues  associated  with  rotating  stall  are  given  by 


rn  = 


2m, , 


2N  +  1 


[1 


respectively,  for  1  <  n  <  N,  where  un  is  an  arbitrary  nonzero  column 
size  2,  Moreover,  tnrn  =  1  if  and  only  if\\un\\  =  ^ujun  =  1/V2- 


0], 

(38) 
vector  of 


Proof.  See  the  Appendix. 


□ 


3ActuaIIy  the  a  parameter  in  [2]  is  different  from  what  we  used  here.  Taking  fi  -*■  tt/(2o)  will 
then  give  exactly  the  same  expression  as  in  [2], 
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With  rn  as  in  Lemma  3.1,  ( rn  ■  r„  =  r„2),  (r„  •  r„)  =  \rn\2,  and  (r„  ■  r„  •  r„)  are 
given  by 


•2  Pn  f 


©2  n 
0 


1 

21V  +  1 


(  M«2  \ 


knl2  =  — 


<?« 1 

_  -  Pn 

L  o  j 

>  rn  •  rn  ‘  rn  —  «— 

2V2  | 

_  0  _ 

,  S„  =  tan-1  —  ,  an 
\Mnl/ 


_  Pnr" 
2  ' 


(39) 

(40) 


To  determine  the  stability  of  the  nth  pair  of  the  critical  modes  corresponding  to 
rotating  stall,  column  vectors  p,n  and  vn  as  defined  in  (75)  must  be  computed.  It  is 
noted  that 


SQo[rn,rn]  =  Y 


,2  r  at 


2B^(2N  + 1  )lc 


'2N+1  t* 


£,J/2r  o 

^(|-i)rT£-'r  o2W+1 


oT 

U2N+1 


'44', 


e#  1 
0 


3//c3/j; 


2W2 


2W2 


« ( h.  _  a  r  °2Af+i  i  g 

Vw  ~j[e-1/2t  \n 


fin  — 


3HC1PZ  2N  + 

4w2  y  ic 


/2AI  +  1 

'  0  ' 

/  zc 

In'  V 

.eN . 

with  y 

=  Yc, 

1 

£1^ 

1  —  1 

lj  S-'D-' 

'  0  ‘ 
ew  _ 

(41) 


37/C3P2  j2N 
AW2 

3Hc3p2 


W  /L 


02W+i  TrE~l/2eN 


l 


?(M[ 


02/v+i  ejv 
2^  0 


2BJ(2N+\)l, 
-1 


A 


AW2lc 

Using  the  expression  for  Do,  we  obtain 

Do 

It  follows  that 

fin 


A 


o 


ft  -  1 

~ 1  1 

R 

0 

1 

* 

\_—k  a/lc_ 

k2  4-  ax /lc 

.  *  *  . 

3  Hcjpl 


xeu 

U  v 

—K 

L  IB  . 

a  —  0. 


a  =  0. 


.ifo] 

0  UJ 

(42) 

(43) 

(44) 


4  W2(k2Ic  +  at) 

The  computation  of  the  vector  v„,  which  is  more  complicated,  is  summarized 
in  the  next  result. 
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Lemma  3.2.  The  vector  vn  as  defined  in  (75)  for' the  local  multimode  model  is 
given  by 

v  _  3Hcipl(c02n  +  2(On  +  Oi2nj)jejUn  _  A  f  ©2 n  1 

4m2nW2(a)ln-4col+ajn-2a2na)nj)\W  )[  0  j’ 

where  1  <  n  <  N /2,  and  ct2n  =  (a  —  ji(2n)2)  •  For  N  <  2n  <  2N,  vn  has 

the  same  expression  as  in  (45)  except  that  2 n  is  replaced  by  t]  =  2 (N  —n)+  1, 
and  ©2 «  fry  ©-,,. 

Proof.  See  the  Appendix.  □ 

The  next  result  gives  the  stability  condition  for  each  projected  dynamics  along 
the  nth  pair  of  the  critical  eigenvectors. 


Theorem  3.3.  For  0  <  2 n  <  N,  the  stability  characteristic  value  (SCV)  for  the 
nth  pair  of  critical  modes  corresponding  to  rotating  stall  is  given  by 


jM  _  3#c3 pi  \W  (&  _  \2  _  WtC3  \ 

2  mnW‘ *  4  \  W  /  \  2  k2Ic+ccz) 


V2n<-o\n  +  Of  |n  +  2a2nC»n^2n 
(co%n~4CO2-ha2n)2  +  4a!2nC02’ 


(46) 

(47) 


w/zere  oi„  =  *.  For  N  <2n  <  2 N.  X^  w  fAe  ifl/ne  as  in  (46)  except  that  2 n 

is  replaced  by  rj  =  2 (N  —  n)  +  1. 
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Re[4(r„  •  vn)]  = 


3Q%n  Hc3P2n 
4m2nW2 


(53) 


Qa'n  =  Im  (^2«  +  2  (On  +  0L2nj)/m2n 
|_<4,  -  4a)n  +  a2n  ~  2a2nCO„j 

Hence  the  expression  for  X^  in  (46)  can  be  verified  following  (76). 


(54) 

□ 


As  mentioned  earlier,  X^  <  0  for  n  =  1,2, N  may  not  imply  stability  of 
the  rotating  stall  dynamics  because  of  the  coupling  between  different  pairs  of  the 
critical  modes.  However  for  the  special  case  of  p  =  0  (that  is,  when  N  pairs  of 
the  critical  modes  in  (30)  cross  the  imaginary  axis  simultaneously),  Xf1’  <  0  for 
n  =  1 , 2, . . . ,  N  determines  the  stability  of  the  rotating  stall  dynamics,  which  is 
surprising.  The  next  lemma  gives  the  expressions  of  p„,  and  vn,  as  in  (79),  which 
will  be  useful  in  establishing  the  final  result  of  our  paper.  The  derivations  of  pn, 
and  vni  are  similar  to  those  of  pn  and  v„,  so  the  proof  is  omitted. 


Lemma  3.4.  Suppose  that  ji  —  0.  Then  the  vectors  pni  and  vni  as  defined  in  (79) 
for  the  local  compression  model  are  given  by 


dm  — 


vni 


3Hc3pnpije^Sn  Sl)  / 

*Pc 

A 

'0 

n—l 

4m\n-i\W2(o)\n-i\  -con  +(oi )  ' 

w 

V 

- 

0 

3Hc3pnpijej(Sn+Sl'> 

f*c 

A 

®tl+l 

4mn+iW2(con+i  -a)„~  (oi) 

7 

L 

0  J 

3Hc3  pn  pi  j  ej  (5" +s,)  /  <pc 

A 

fe- 

-ul 

4 mvW2(a)q  -con~wi)  \W 

7 

L  0  . 

’ 

(55) 


0  <  n  + 1  <  N, 

N  <  n  +  i  <  IN , 
(56) 


where  n,  l  =  1,  2, . . . ,  N,  n  ±  l,  r)  =  2 N  -f  1  —  n  —  l,  <pc  is  as  in  (36),  and  0„  is 
as  in  (37). 


Theorem  3.5.  Under  the  conditions  p  =  0,  c3  <0,  and  ci  +  3c3  =  0,  the 
rotating  stall  dynamics  are  locally  stable  if  and  only  if  cq  +  10c  3  >  0. 


Proof.  Local  stability  of  rotating  stall  implies  that  each  pair  of  the  corresponding 
critical  modes  is  locally  stable,  which  implies  that  X^  <  0  forn  =  1,2 , . . . ,  N. 
Because  <pc  =  2 W,  p  =  0,  and  a  =  0  at  the  criticality,  a2n  =  0  for  1  <  n  <  N.  It 
follows  that  for  n  =  1, 2, . . . ,  N,  k(2n)  <  0  is  equivalent  to 


3^3#;  (l_  3Hxc3\ 
W3  \4  Wlctc2  ) 


(57) 


Because  C3  <  0  and  all  other  parameters  are  positive,  this  relation  is  equivalent  to 


1  3Hxc3 

- _  >  0 

4  WIck2 


(58) 
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Noting  that  the  expression  of  r  at  criticality  becomes 

W 

~4  B2H!c(c0  +  c1+c3)’  (59) 

straightforward  calculation  proves  that  X^  <  0  is  equivalent  to  co  +  IOC3  >  0. 

For  sufficiency,  assume  that  co  +  IOC3  >  0,  which  is  equivalent  to  X^  <  0  for 
n  =  1, 2, . . . ,  N.  It  will  be  shown  using  Theorem  5.3  that  the  compressor  model 
is  locally  stable  at  uniform  flow  equilibrium  <pe  =  cpc,  and  \jre  =  t/rc  with  y  =  yc. 
Indeed  there  holds  Xnn  =  SA.^  <  0  for  n  =  1,2 , . . . ,  N.  Furthermore,  through 
straightforward  calculations, 


4 (h  ■  vni)  =  jrvmnPnPl  (2 N  +  1), 

O 

4 (Xi  ■  Hm)  =  gn-~(2N  +  1), 

for  n  ^  l  where  r;x,  rv  are  real  numbers.  Thus 

Re  [4 •  Mn;)]  =  Re  [4 (if  •  vn/)]  =  0 
for  n,  l  ==  1,  2, . . . ,  N.  It  follows  that 


Xni  =  16Re  Q[rn,pt]  +  j€„C[r„,r/,  r/]^ 
at  y  =  yc.  Using  the  expression  of  r„,  and  (pc  =  2 W, 


rn-  Ml  = 


3Hc3xpj 
"  4W2k2Ic 
By  (48),  (49),  and  relation  £nrn  =  1, 


rn. 


n-n 


pf 


rn- 


(60) 

(61) 

(62) 

(63) 

(64) 


4e[rn,/zz]  =  -zVf^7^")  -  4C[r„.rrF/]  =  |^£.  (65> 


mnlc  \  2W2k  ) 

The  expression  of  Xni  is  then  obtained  as 

12 Hc3pf  3Hxc3\ 


Xni  = 


m„  W3 


i 1  3Hxc3\ 

\2  WICK2)' 


(66) 


Again,  because  C3  <  0  and  all  other  parameters  are  positive,  Xni  <  0  is  equivalent 
to 

1  3Hxc3 

<67> 

which  is  implied  by  (58).  Hence  the  bifurcated  system  is  locally  stable  in  light  of 
Theorem  5.3.  □ 


It  is  very  surprising  to  see  that  the  stability  condition  for  rotating  stall  is  identical 
to  the  case  of  N  =  1  for  the  third-order  Moore-Greitzer  model  [9],  However,  it 
should  be  noted  that  the  result  of  Theorem  5.3  is  derived  based  on  assumptions 
that  p  =  0,  and  x/rc(  )  is  cubic  satisfying  C3  <  0  and  ci  +  3c3  =  0  as  in  [10],  [9]. 
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4.  Conclusion 

Rotating  stall  dynamics  were  investigated  in  this  paper  for  axial  flow  compres¬ 
sors,  and  necessary  and  sufficient  conditions  were  obtained  for  local  stability  of  the 
rotating  stall  dynamics  in  terms  of  various  compressor  parameters.  Although  our 
stability  conditions  were  derived  for  the  multimode  Moore-Greitzer  model,  they 
hold  for  the  limiting  case,  i.e.,  the  full  PDE  Moore-Greitzer  model.  This  result  is 
important  because  past  research  has  almost  exclusively  focused  on  the  third-order 
Moore-Greitzer  model  ( N  =  1),  of  which  many  active  control  laws  were  proved 
to  be  effective  in  the  literature  but  fail  to  work  in  engineering  practice  because 
of  the  spillover  of  the  high-order  modes  [4],  Our  local  stability  conditions  clearly 
offer  new  insight  into  the  design  and  active  control  of  rotating  stall  for  the  mul¬ 
timode  Moore-Greitzer  model,  thereby  improving  the  performance  of  axial  flow 
compressors. 


5.  Appendix 

This  section  reviews  briefly  the  projection  and  Lyapunov  methods  for  local 
stability  analysis  of  Hopf  bifurcations  and  gives  some  of  the  proofs  for  Section  3. 

For  stability  analysis  of  Hopf  bifurcations,  the  system  under  consideration  is 
the  following  nth-order  parameterized  nonlinear  system: 

x  —  f(y,  x),  /(y,0)  =0  Vy  e  R  and  x0  €  5  C  R",  (68) 

where  x  e  R",  y  is  a  real- valued  parameter,  and  S  is  a  linear  subspace  of  R", 
to  be  clarified  later.  It  is  assumed  that  /(•,  •)  is  sufficiently  smooth  such  that  the 
equilibrium  solution  xe,  satisfying  / (y,  xe)  =  0,  is  a  smooth  function  of  y.  The 
smoothness  of  /(■,  •)  implies  the  existence  of  a  Taylor  series  expansion  near  the 
origin  of  R"  in  the  form 

x  =  f(y,x)  =  L(y)x  +  Q(y)[x,x ]  +  C(y)[x,x,x]  -1 - ,  (69) 

where  L(y)x,  Q(y)[x,  x],  and  C(y)[x,  x,  x]  can  be  expanded  into 

L(y)x  =  Lox  +  8yL\x  +  8y2L2X - ,  (70) 

Q(y)[x,x]  =  Qo[x,x]  +  SyQi[x,x]  +  ---,  (71) 

C(y)[x,  x,  x]  =  Co[x,  x,  x]  +  SyCilx,  x,  x]  -i - ,  (72) 

with  Sy  =  y  —  yc,  and  Lo,  L\,  Lj  constant  matrices  of  size  (n  x  n).  Suppose  that 
L{y)  possesses  m  (<  n/2)  pairs  of  complex  eigenvalues  Xkiy)  —  ajt(y)db  jpkiy), 
dependent  smoothly  on  y.  It  is  assumed  that  for  1  <  k  <  m, 

,  da* 

otkiYc)  -  0,  P(yc)  =cok^0,  cck(yc )  =  (Yc)  ±  0,  (73) 

dy 

while  all  other  eigenvalues  of  L(y)  are  stable  at  and  in  a  neighborhood  of  y  =  yc . 
That  is,  m  pairs  of  eigenvalues  cross  the  imaginary  axis  simultaneously.  Then  each 
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A*(y)  is  a  critical  eigenvalue,  and  so  is  its  conjugate.  The  center  space,  or  the 
eigenspace,  for  the  m  pairs  of  critical  eigenvalues  is  denoted  by  Sc,  and  is  assumed 
to  be  orthogonally  complementary  to  S  in  the  sense  that  S  ©  Sc  =  R" .  Because  the 
projection  of  *0  €  S  to  Sc  is  zero,  the  equilibrium  xq  satisfying  (68)  will  be  called 
the  zero  solution.  The  strict  crossing  assumption  implies  that  the  zero  solution 
changes  its  stability  as  y  crosses  yc.  For  instance,  a'k(yc)  >  0  implies  that  the 
zero  solution  is  locally  stable  for  y  <  yc  and  becomes  unstable  for  y  >  yc.  The 
crucial  problem  is  the  determination  of  local  stability  near  the  critical  parameter 
yc  at  which  Hopf  bifurcations  occur  and  the  periodic  solutions  are  born. 


A.  1.  The  projection  method. 

The  center  space  Sc,  spanned  by  all  the  critical  eigenvectors,  is  complete.  If  all  the 
critical  eigenvalues  are  distinct,  then 

Sc  =  SC1  ©  SCl  ©  ■  •  ■  ©  SCm ,  (74) 

where  SCk  denotes  the  eigenspace  spanned  by  the  Ath  pair  of  the  critical  eigen¬ 
vectors.  The  projection  method  in  [7]  projects  the  nonlinear  dynamics  into  the 
subspace  of  SCk  so  that  stability  of  the  projected  dynamics  can  be  analyzed.  Asso¬ 
ciated  with  each  pair  of  the  critical  eigenvalues  (A*,  A*),  there  can  define  an  SCV 
(stability  characteristic  value)  A®  whose  sign  determines  local  stability  of  the 
projected  dynamics.  See  [7],  [6],  [1]  for  details.  Let  the  Taylor  series  of  f(y,  x) 
be  of  the  form  in  (69),  where  Lo  =  L( 0).  An  algorithm  to  compute  A®  is  now 
outlined. 

Step  1 .  Compute  the  left  row  eigenvector  l k  and  the  right  column  eigenvec¬ 
tor  r/c  of  Lo  corresponding  to  the  kth  critical  eigenvalue  of  A*(0)  =  jwk. 
Normalize  by  setting  lkrk  =  1. 

Step  2.  Solve  column  vectors  fik  and  vk  from  the  equations 

-L0Pk  =  ~  Qo[rk,  OtL  (2jcokI  -  L0)vk  =  rk\.  (75) 

Step  3.  The  coefficient  A®  is  given  by 

A®  =  2Re  {240o[r*,  /**]  +  tkQolh,  «*]  +  ^4C0[r*,  rk,  rk]  j  .  (76) 


Theorem  5.1.  Suppose  that  Lo  has  N  pairs  of  nonzero  critical  eigenvalues  on 
the  imaginary  axis  with  the  rest  on  the  open  left  half-plane.  Then  the  projected 
dynamics  along  the  kth  pair  of  the  eigenvectors  is  stable  ifl.2  <  0,  and  unstable 
if  A 2  >  0.  For  the  case  ofN  =  1,  the  associated  Hopf  bifurcation  is  supercritical 
or  stable  ifk 2  <  0,  and  subcritical  or  unstable  ifki  >  0. 
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It  should  be  clear  that  if  Lo  admits  only  a  single  pair  of  imaginary  eigenvalues, 
then  local  stability  of  the  bifurcated  system  is  equivalent  to  X®  <  0.  However, 
if  Lo  admits  more  than  one  pair  of  imaginary  eigenvalues,  local  stability  of  each 
projected  dynamics  does  not  guarantee  local  stability  of  the  bifurcated  systems.  In 
this  case,  local  stability  of  the  bifurcated  systems  can  be  tackled  by  the  Lyapunov 
method  [3]. 


A.  2.  The  Lyapunov  method. 

The  Lyapunov  method  for  multiple  critical  modes  has  been  developed  by  Fu  [3]. 
The  method  gives  sufficient  conditions  on  the  existence  of  the  Lyapunov  function  to 
guarantee  local  stability  of  the  bifurcated  system.  The  following  notion  is  needed. 

Definition  5.2.  Under  the  condition  that  Lo  admits  N  >  1  pairs  of  eigenvalues  on 
the  imaginary  axis  with  the  rest  on  the  open  left  half-plane,  the  nonlinear  system 
(69)  is  said  to  be  locally  center-symmetric  in  the  sense  of  Lyapunov  if  the  following 
conditions  hold: 

coj  ^  a >m,  cot  ^  2wm  cot  ^  3 com  whenever  N  >  2, 

ojj  /  com  +  a>k,  cot  2a>m  -j-  cok,  2m,-  ^  com  +  <ok,  whenever  N  >3, 

a>i  com  +  (Ok  +  m,  ^  +  com  ^  cok  +  coi,  whenever  N  >  4, 

hold  where  all  indices  are  distinct. 

Theorem  5.3.  Suppose  that  Lo  admits  N  >  1  pairs  of  eigenvalues  on  the  imagi¬ 
nary  axis  with  the  rest  on  the  open  left  half-plane,  and  the  nonlinear  system  (69)  is 
locally  center-symmetric  in  the  sense  of  Lyapunov.  Then  there  exists  a  Lyapunov 
function  that  guarantees  local  stability  of  (69)  at  the  criticality  if 

Xkk  =  16Re  4  ^2Q0[rk,  MkJ  +  QolLk,  He]  +  ^C0[rk,  ik,  rk]^  J  <  0, 

for  k  =  1 ,  •  •  • ,  N,  (77) 

Xkp—  16Re  1 4  ^Qo[rk,  /xp]  +  Q0jrp,  i^]  +  Q0[rp,  Mk|J  +  ^Q)[rk,  rp,  rp|j  j  <  0, 

(78) 

for  k,  p  =  1, . . . ,  N  with  k  p,  where  vk  are  the  solutions  to  (75),  and  pkp, 
vkp  are  given  by 


1 

Pkp  =  --(L0  -  j(cok  -  (Op))  Qoln,  rp], 
Vkp  =  ~^(Lo  -  j(cok  +  (op)rl  Qdrk,  rp ]. 


(79) 
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It  is  noted  that  Xkk  =  8k®-  Hence  local  stability  of  Hopf  bifurcations  with 
multiple  pairs  of  critical  modes  requires  more  than  local  stability  of  each  pro¬ 
jected  dynamics.  Roughly  speaking,  the  nonpositivity  of  Xkp  with  k  ^  p  takes  the 
coupling  of  different  pairs  of  critical  modes  into  account  for  local  stability  of  the 
bifurcated  system  in  (69). 

A3.  Proof  of  Lemma  3.1. 

Note  that  D  =  SLS _1  as  in  (27).  If  ro  is  an  eigenvector  of  D,  then  S~l  ro  is  an 
eigenvector  of  L.  Let  the  nth  eigenvector  of  D  be  denoted  by  xn  +  jyn  where  xn 
and  yn  are  real  vectors  of  size  2 (N  +  1)  corresponding  to  the  eigenvalue  of  ja>n 
at  y  =  Yc  where  con  —  n/(bmn)  >  0  for  n  =  1, 2, . . . ,  N.  By  equating  real  and 
imaginary  parts  of  the  eigenequation  for  D,  the  two  equations 

Dxn  =  -wnyn,  Dy„  =  wnxn,  (80) 

are  obtained.  This  yields  the  eigenequations 

D2xn  =  —w2x„,  D2yn  =  —(o2yn.  (81) 

Hence  xn  and  yn  satisfy  the  same  eigenequation.  We  need  to  solve  only  one  of 
them;  the  other  can  be  obtained  from  either  Dxn  =  —conyn  or  Dyn  —  conxn.  At 

y  =  Yc, 

D2  =  diagCC^,  D2, . . . ,  D2n),  D2  =  —co2I2,  n  =  \,2,...,N.  (82) 

Thus  it  follows  that  the  (2  x  2)  matrix  Dn  at  the  (n  +  l)th  diagonal  block  of  D2 
is  canceled  in  the  sum  of  D2  +  w2I.  Because  (D2  —  co2I)xn  =  0,  and  det(Dg)  = 
(k2  +  xa/lc )2  ^  0,  xn  has  the  form 

Xn  =  [02n  02(N-n)  ]T .  ®2k  €  ^ ,  (83) 

where  un  is  an  arbitrary  nonzero  real  column  vector  of  size  2.  Now  yn  can  be 
calculated  as 

l  1 

y«  —  Dxn  —  [  O2 n  (DnUn)^  O2 (N—n)  ]  =  [02/i  vj  ^2(N—n)  ]  • 

a>n  con 

(84) 

The  right  eigenvector  rn  corresponding  to  the  eigenvalue  jcon  is  thus  given  by 


r  O2  n 

r, I  =  s~\xn  +  jyn )  =  S-1  un  +  jv„ 

r  0 

’  Vn  ~  L-1 

r 

of- 

(85) 

-  ®2{N-ri) 

for  n  =  1,  2, . . . ,  N.  The  relation  between  un 

and  v„  implies  that 

un+jv„=jun,  j=  y 

7]. 

(86) 
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where  n  =  1, 2, . . . ,  N.  It  is  now  clear  that 


0 2N+1 

tte-i/2- 

o2n 

_  {2BJ{2N  +  l)lc)~l 

0T 

u2Af+l  J 

Un  +  JVn 
-02(JV-»-t)  - 

I  2mn  1 
2N  +  1 


cosn#i 


sin  n6\ 


cos  n02N+i  sin  n@(2N+i) 
0  0 


(Un  + 


(87) 


which  becomes  the  same  as  (38)  by  substituting  the  expression  of  un  +  jvn  as  in 
(86).  Next  denote  Cd  as  the  left  eigenvector  of  D.  Because  of  the  block  structure 
of  D  and  Dn  =  — Z)J  at  y  =  yc  for  n  =  1, 2, . . . ,  N,  Id  —  r^.  Hence 

ln  ~  ^  ~rD*-rn>  *  ~  4  [  0Tw+i  4B2lc(2N  +  !)  _ 
which  gives  the  expression  of  ln.  The  normalization  condition  tnrn  =  1  yields 
4  rn  =  mnr^rn  =  2||h„||2  =  1,  (89) 

which  is  equivalent  to  ||  «„  H  —  1  /\/2.  The  proof  is  now  complete. 


=  r„m„,  (88) 


A.4.  Proof  of  Lemma  3.2. 


By  definition  in  (75), 

vn  =  ~  Lo)~'  Qolrn ,  rn]  =  -  Dy'SQolrn,  rn] 


2 

=  ^S-\ljwn-D)-lS 


)rTzr1r  o2W+1 


oT 

U2W+1 


44’r  J 


®2«  L;2«„ 

0  j  ‘ 

(90) 


There  are  two  cases  to  consider: 


(i)  2 n<N 

(ii)  N  <2n<  2  N. 

For  (i), 


©2  n 

0 

3 Hc3pleJ2Sn  (<pc  \ 

'  cos(2n0i)  s 

iin(2n6*i) 

:  1 

4m2nVk2  \W  J 

cos(2ft<92iV+i)  sin(2w02A’+i)  j 

0  0  J 
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x  (2 jwn  -  Dln)  1 
At  the  critical  parameter  y  =  yc, 

n  —  0  -con  1 

”  L  con  0  J 

Then  straightforward  computation  gives 


u>n 


-i 
n  ' 


(91) 


(92) 


(2 j(On  -  Dln)  1 


1 

(C02n  +  2  COn  +  CC2  nj)j 

1 

(p In  -  4 oil  +  a\n  -  2a 2n(onj 

(93) 


Hence  for  n  <  N / 2,  direct  computation  gives  the  expression  in  (45). 

For  (ii),  without  loss  of  generality,  it  is  assumed  that  6k  =  2 (k  —  1)tt/(2 N  + 1), 
where  1  <  k  <  2 N  +  1.  Because  N  <  2n  <  2 N, 


2 n  =  2N  +  1  —  rj,  0  <  r]  <  N. 

It  follows  that  for  1  <  k  <  2N  +  1, 

e~j2n8k  _  e~j(2N+l-ri)6k  _  £ji]8k , 

where  0  <  ij  <  N .  In  this  case, 


(94) 

(95) 


vn  = 

Noticing  that 

(2;ai„  -  D„)_1  - 


(2jwn  -  D„)  1 


3 Hc3pleJ2Sn  (4>t 


4m,,  W2 


(f-1) 


5-1(2>„-D)-15 


-l , 


-  4a)n  +an~  2aVMnj 


2  jcon  -  an 


0 


(96) 


-co„ 


(Or, 


2  jcon  -  an 


V 

(&>„  +  2oj„  +  a„j)j 

V 

J . 

~  4 a>l  +afj-  2anconj 

J  _ 

(97) 

(98) 


direction  computation  gives  the  same  vn  as  in  (45)  except  that  2 n  is  replaced  by 
—t]  with  a)-,,  =  co,,,  a -n  =  an,  and  m -n  =  m 


References 

[1]  E.  H.  Abed  and  Jyun-Homg  Fu,  Local  feedback  stabilization  and  bifurcation  control,  I  —  Hopf 
bifurcation.  Systems  Control  Lett.,  7,  11-17,  1986. 

[2]  R.  A.  Adomaitis  and  E.  H.  Abed,  Bifurcation  analysis  of  nonuniform  flow  patterns  in  axial-flow 
gas  compressors,  1st  World  Congress  of  Nonlinear  Analysis,  Aug.  1992. 

[3]  J.-H.  Fu,  Lyapunov  functions  and  stability  criteria  for  nonlinear  systems  with  multiple  critical 
modes,  Math.  Control  Signals  Systems,  1,  255-278,  1994. 

[4]  G.  Gu,  S.  Banda,  and  A.  Sparks,  An  overview  of  rotating  stall  and  surge  control  for  axial  flow 
compressors,  Proc.  35th  IEEE  Conf.  on  Dec.  and  Contr.,  Kobe,  Japan,  2786-2791,  Dec.  1996. 


350  Gu,  Sparks,  and  Belta 


[5]  R.  A.  Horn  and  C.  R.  Johnson,  Topics  in  Matrix  Analysis,  Cambridge  University  Press,  London 
and  New  York,  1991. 

[6]  L.  N.  Howard,  Nonlinear  oscillations,  in  F.  C.  Hoppensteadt,  ed..  Nonlinear  Oscillations  in 
Biology,  American  Mathematical  Society,  Providence,  RI,  1979,  pp.  1-68. 

[7]  G.  Iooss  and  D.  D.  Joseph,  Elementary  Stability  and  Bifurcation  Theory ,  Springer- Verlag,  Berlin 
and  New  York,  1980. 

[8]  C.  A.  Mansoux,  J.  D.  Setiawan,  D.  L.  Gysling,  and  J.  D.  Paduano,  Distributed  nonlinear  modeling 
and  stability  analysis  of  axial  compressor  stall  and  surge,  Proc.  1994  Amer.  Contr.  Conf,,  2305- 
2316, 1994. 

[9]  F.  E.  McCaughan,  Bifurcation  analysis  of  axial  flow  compressor  stability,  SIAM  J,  Appt.  Math., 
vol.  20, 1232-1253,  1990. 

[10]  F.  K.  Moore  and  E.  M.  Greitzer,  A  theory  of  post-stall  transients  in  axial  compressors,  Part  1  — 
Development  of  the  equations,  ASME  J.  of  Engrg.  Gas  Turbines  and  Power,  vol.  108,  68-76, 
1986. 


